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Complex quantum systems consisting of large numbers of strongly coupled states exhibit char¬ 
acteristic level repulsion, leading to a non-Poisson spacing distribution which can be described by 
Random Matrix Theory. Scattering resonances observed in ultracold atomic and molecular systems 
exhibit similar features as a consequence of their energy level structure. We study how the over¬ 
lap between Feshbach resonances affects the distribution of resonance spacings. The spectrum of 
strongly overlapping resonances turns out to be non-Poisson even when the assumptions of Random 
Matrix Theory are not fulhlled, but the spectrum is also not completely chaotic and tends towards 
being semi-Poisson. 
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Introduction. Magnetic Feshbach resonances are 
widely used to control atomic interactions in ultracold 
quantum gases [T]. They have been extensively studied 
theoretically and experimentally for many different com¬ 
binations of alkali-metal-atoms, for which detailed theo¬ 
retical models have been constructed to understand the 
Feshbach spectra that result from tuning the magnetic 
field. In these cases quantum numbers can readily be 
assigned to the molecular states that give rise to indi¬ 
vidual resonances. Precise experimental control and the¬ 
oretical understanding combine to allow the production 
of weakly bound Feshbach molecules [2]. The STIRAP 
technique [3] , which allows adiabatic transfer of the pop¬ 
ulation to deeply bound states by using properly designed 
laser pulses, has opened the possibility to study ultracold 
molecular gases [4HT0]. 

Close to an isolated Feshbach resonance caused by a 
ramping bound state that crosses the collision thresh¬ 
old at energy E = 0 as magnetic field B is tuned, the 
s-wave scattering length for a near-zero-energy collision 
that has no inelastic loss is described by the well-known 

formula mm 

a{B) = abg(l- ^ y (1) 

where A is the resonance width and is the position 
of the resonant pole in a{B). However, it is often the case 
that two or more Feshbach resonances overlap with each 
other so that a proper description of the scattering prop¬ 
erties requires taking all the relevant bound states into 
account mm- This is most pronounced for non-alkali 
systems such as highly magnetic lanthanide atoms, where 
anisotropic interactions introduce strong mixing of many 
partial waves, which can result in hundreds of overlap¬ 
ping Feshbach resonances |15[ I16j . In experiments per¬ 
formed with erbium m and dysprosium nails] atoms, 
the observed resonance density is of the order of a few 
resonances per gauss for bosonic atoms with no nuclear 
spin and exceeds 20 resonances per gauss for the fermion 
^®^Er, which has additional hyperfine structure. In such 


cases there is no way to assign quantum numbers to in¬ 
dividual resonances or to predict their positions theoret¬ 
ically. Interestingly, the resonances can still be used to 
produce weakly bound dimers [20] . 

Statistical analysis of the Feshbach spectra recorded 
for erbium m showed that the spacings between 
neighboring resonances exhibit characteristic correla¬ 
tions. While one would expect a Poisson distribution of 
spacings for uncorrelated energy levels, the experimen¬ 
tal results were closer to a Wigner-Dyson distribution 
(WD) [3TJ[52|, which indicates strong repulsion between 
the overlapping resonances. This kind of level spacing 
distribution can be derived by means of random matrix 
theory (RMT) and has been shown to accurately repro¬ 
duce resonances in nuclear reactions [531IHI • Recently, 
it has been predicted that possible resonance densities in 
ultracold atom-molecule and molecule-molecule collisions 
can be so high that they may exhibit similar behavior. 
References [53 [53 show how to estimate the density of 
states and use RMT to construct a statistical scattering 
model. The Bohigas-Giannoni-Schmit conjecture implies 
that for systems with such statistical properties the clas¬ 
sical dynamics would be chaotic m- 

On the other hand, the scattering length for an arbi¬ 
trary set of overlapping Feshbach resonances can also be 
found analytically, e.g. by using using methods of mul¬ 
tichannel quantum defect theory (MQDT) [T3]. In this 
paper, we apply this analytical model to the case of ex¬ 
tremely dense Feshbach spectra. Interestingly, we find 
that the distribution of the resonance spacings starts to 
resemble the Wigner-Dyson distribution when the mean 
resonance width exceeds average spacing between the res¬ 
onances, even if the RMT assumptions are not fulfilled. 

Scattering length in the presence of many resonances. 
We start our analysis by considering a general scattering 
problem of two particles with complex internal structure 
giving rise to a large possible number of scattering res¬ 
onances. Instead of studying in detail the origin of the 
resonances as done in |15| . we will simply assume that 
there are iV ^ 1 bound states of the composite 2-particle 
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system which can cross with the single entrance channel 
as B is tuned over a range spanned by these levels (it is 
also possible to consider fields other than a magnetic one 
to do the tuning). This channel is characterized by its 
background scattering length atg- We will assume that 
the particles interact via a van der Waals potential at 
long range and express atg as ra with a being the mean 
scattering length of the van der Waals potential [55j . The 
characteristic range of interparticle separation at which 
interchannel couplings are present will be assumed to be 
much smaller than a. Under these assumptions, the mag¬ 
netic field dependence of the scattering length for a set of 
N overlapping Feshbach resonances can be conveniently 
described by the following analytical formula [13] 


a{B) =abg [l- 


N 

E 


A,. 


■^^b-b,-sb,-y: 


B — Bj X D 


( 2 ) 

Here is the local width associated with fc-th reso¬ 
nance, Bk is the position at which the bare fc-th bound 
state crosses the threshold and 6Bk is the corresponding 
shift of the resonance position due to interaction with 
the open channel, for van der Waals interactions given by 
This formula can be derived using 
the MQDT formalism described in detail in m, where 
it was able to account for overlapping resonances in sev¬ 
eral different examples involving alkali-metal atoms. The 
actual positions of resonances can be calculated by 
finding the zeros of the denominators in Eq. In gen¬ 
eral, the widest and most closely lying resonances tend 
to have the strongest impact on the other ones. 

In the derivation of this formula, the closed channels 
are chosen to be uncoupled, with the only nonzero cou¬ 
plings existing between the entrance channel and the 
closed channels. Another assumption of our model is 
linear variation of the tangent of the bound state phase 
ta,ni/k{E) oc E—Sfi{B—Bk) for each closed channel state. 
As a result of these approximations, our model is appli¬ 
cable only near the energy threshold. Formula ([^ de¬ 
scribes the scattering length in terms of “local” resonance 
widths. The fc-th width A^ is determined only by the 
coupling with the bound state crossing the threshold at 
Bk- It is also possible to algebraically transform Eq. (§ 
to the form 


= (3) 

Here the resonance widths are a function of all interchan¬ 
nel couplings instead of local parameters. Regardless of 
the choice, the resonance positions are the same and it is 
the distribution of B^es and the spacings between them 
that is of our current interest. 

Dense sets of resonanees. Our MQDT formalism suc¬ 
cessfully reproduced the complex structure of a{B) when 
several overlapping resonances occur in the case of alkali 
metal species m- It can also be used to set up a test sys¬ 
tem for a dense set of resonances in which we can examine 


the change in level statistics that occurs as the coupling 
strength increases, with different assumptions about the 
distribution of the interactions A^, or of the “bare” cross¬ 
ing positions Bk before coupling to the entrance channel 
is turned on. Let us emphasize that we are not intending 
to produce realistic models for atomic species like Er or 
Dy, since our model of the potentials and interactions is 
too simple to capture the physics of such complex atoms. 
Rather, we are interested in understanding some of the 
conditions that can give rise to Poisson or WD distribu¬ 
tions when we have a set of interacting tunable Feshbach 
resonances in cold atomic collisions. 

The conventional way of dealing with systems with 
large number of coupled channels with no good quantum 
numbers (apart from time-reversal symmetry) is to use 
RMT. In this approach all the couplings are assumed to 
be random numbers with statistical properties governed 
by Gaussian Orthogonal Ensemble. As a result, the en¬ 
ergy level spacings as well as resonance spacings follow 
the Wigner-Dyson distribution, defined as 

Pm/d(s) = yexp , (4) 

where pwd{s) gives the probability density and s is the 
distance between adjacent levels in units of mean reso¬ 
nance spacing d = 1/p with p = N/Bmax, where we as¬ 
sume B to vary from 0 to B^ax for convenience. While 
normally the Hamiltonian of the system is fixed and the 
distribution of energy levels is characterized, in studying 
Feshbach resonance spacings we fix A = 0 and look at 
the resonance distribution as the Hamiltonian is varied 
by changing B. 

Within our treatment we start from Eq. and make 
different statistical assumptions about the bare resonance 
positions and coupling strengths. In the simplest scenario 
one can assume all couplings to be equal and all bare 
resonance positions to be completely uncorrelated. The 
distribution of spacings between adjacent bare levels is 
then a Poisson one that can be described by 

Pp{s) = e"* (5) 

This clearly violates the assumptions of RMT. The cor¬ 
relations between observed resonance positions will come 
solely from the terms SBj in Eq. (H) that describe 
mutual interaction between the resonances. For very 
weak couplings, resulting in mean A (and consequently 
SB) much smaller than average distance between the res¬ 
onances, the distribution of spacings should intuitively 
stay a Poisson one. However, for large widths the reso¬ 
nances will start to overlap significantly and the statistics 
will be modified. For a given set of resonance parameters 
{r,Ai, Bi} we obtain the statistics of resonance spacings 
numerically and fit it to the Brody distribution, defined 
as [29] 

Pb{s) = a{f] + l)s'’exp , 


( 6 ) 
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FIG. 1. (color online) Brody parameter fitted to a set of 
resonances with width A in units of mean spacing, averaged 
over 25 realizations, along with the standard deviation for 
a single realization (shaded regions). Blue: Poisson initial 
distribution of resonance positions, red: Wigner-Dyson initial 
distribution of resonance positions. 


where ry is a real number 0 < 77 < 1 , Q: = r j 
r is the Gamma function. This is a heuristic distribution 
that in the limit ?7 —f 1 this function approaches WD dis¬ 
tribution, while at 77 0 it recovers Poisson statistics. 

Here we treat 77 as a fitting parameter. Then we aver¬ 
age 77 over many numerical realizations corresponding to 
different sets of Bi. 

Results and discussion. Figure summarizes our hnd- 
ings for an example with 50 resonances of equal widths 
A, showing the dependence of Brody parameter 77 on the 
width of the resonances. For uncorrelated bare resonance 
positions (blue line), we observe that at small widths (less 
than the mean spacing d), the spacing distribution is es¬ 
sentially a Poisson one, but in the limit of large widths 
A d it tends towards a WD distribution with the value 
of r] being around 0.7. It is worth noting that this is close 
to the value measured in experiments with ultracold Er 
and Dy [niEQ], although our model assumptions are not 
designed to represent these species. We checked that this 
feature is independent of abg ■ This approach of 77 towards 
values of around 0.7 still holds if we assume a distribu¬ 
tion of initial resonance widths, although when one starts 
from a set containing many very narrow resonances which 
do not influence each other strongly, the effect will be nat¬ 
urally weaker and will require much larger average A. On 
the other hand, if the bare resonance positions come from 
a WD distribution, the additional effect of the coupling 
between resonances due to the width A in ([^ changes 
the statistical properties only very slightly, as shown by 
the red line in Figure 

It is interesting to study the scattering length a{B) as 
well as finite energy scattering properties for an exam- 


FIG. 2. Scattering length versus B in units of the mean spac¬ 
ing p for an example with 50 resonances having equal widths 
of twice the mean spacing. Not every resonance is visible as 
a divergence of a due to finite resolution of the calculated set 
of points. 



FIG. 3. (color online) sin^ d{E) for the same set of resonances 
as on Fig. assuming random magnetic moment differences 
ranging from 0 to 0.5 Eq per magnetic field unit. 


pie case. Figure shows the scattering length for the 
case of 50 strongly overlapping resonances with an ini¬ 
tial Poisson distribution and equal widths. Note that the 
widths of the interacting resonances are not the same, but 
show considerable variation due to interactions among 
the resonances. When the initial resonance width ex¬ 
ceeds strongly the mean resonance spacing d, the average 
observed resonance width should approach the value of 
d. As a result, many narrow resonances can be seen on 
Figure For weakly overlapping resonances the widths 
are not modified so strongly. In contrast to typical situa¬ 
tions in alkali gases, it is hard to extract a unique value of 
Obg from the figure. In fact, we find that the mean scat¬ 
tering length calculated over a range of magnetic field 
will take a value of the order of d, the mean scattering 
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length of the van der Waals potential. Interestingly, two 
recent experiments with Dy atoms naisi] both reported 
finding background scattering lengths having the order 
of magnitude of this mean scattering length. 

The MQDT theory of Ref. [13] analytically continues 
the bound state spectrum and threshold properties into 
the finite E scattering continuum, where resonances can 
cross one another and the resulting scattering phase shift 
can exhibit a complex pattern. An additional parameter 
6fj, for each resonance, describing the magnetic moment 
difference between the open and closed channel states, is 
needed for this extension. Figure shows the value of 
the sine square of the scattering phase shift, sin^^(£'), 
at different energies and magnetic fields for the same set 
of resonances as on Figure]^ assuming Sfi parameter to 
be a random number between 0 and 0.5 Eq per magnetic 
field unit. The Figure shows that at finite energies the 
phase shift varies rapidly with both energy and magnetic 
field as narrow continuum levels cross one another. Such 
behavior would have notable consequences, making the 
scattering properties of a thermal gas (characterized by 
wide distribution of collision energies) much harder to 
control. 

In addition to the level spacing distribution, the fluc¬ 
tuations in the number of resonances among different 
magnetic field intervals are commonly studied for com¬ 
plex systems. They can be measured by the num¬ 
ber variance S^(Ai?) = N'^{AB) — N{AB) , where 
N{AB) is the number of resonances found in the range 
from Bq to Bq + AB and the averaging is performed 
over Bq. Figure shows the results obtained for our 
model for different resonance widths. The spectra are 
far more rigid than a Poisson distribution would pre¬ 
dict but are still quite far from the predictions of WD 
statistics. In fact, we found that the semi-Poisson dis¬ 
tribution, psp{s) = 4se“^®, for which the number vari¬ 
ance S^(Ai?) = ABj2 -I- (1 — |32|, character¬ 

izes the statistical properties of our model much bet¬ 
ter than the Poisson or Wigner-Dyson one. The semi- 
Poisson distribution describes the statistics directly at 
the metal-insulator phase transition in condensed matter 
systems [351 [33]. Fitting a sample of random numbers 
with semi-Poisson distribution to the Brody distribution 
gives the rj parameter of around 0.6, in qualitative agree¬ 
ment with the behavior we observed on Fig. 

Conclusions. We have extended an analytic MQDT 
model that successfully explains sparse sets of overlap¬ 
ping resonances of alkali-metal species to examine the 
statistics of level spacings of dense sets of many overlap¬ 
ping resonances. While such a model is not designed to 
capture the physics of complex atomic species like Er or 
Dy, it demonstrates the onset of a level spacing distri¬ 
bution similar to the chaotic one. If a non-interacting 
set of resonances described by a Poisson distribution of 
spacings is coupled to an ultracold s-wave entrance chan¬ 
nel with a set of widths that are small compared to the 
mean level spacing, then the interacting set maintains 
a near-Poisson distribution. However, when the widths 
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FIG. 4. (color online) Number variance of the resonance 
positions for five different resonance widths, A = 0.2 (blue), 
0.5 (red), 1 (magenta), 1.5 (purple) and 2.5 (green) in units 
of mean spacing d, averaged over 25 realizations. The dashed 
lines show the limiting analytical results for WD (lower line) 
and Poisson (upper line) statistics. The dotted line depicts 
the semi-Poisson distribution. 


begin to exceed the mean level spacing, the level spac¬ 
ing distribution changes, with a Brody parameter rj that 
has a value of around 0.7 and is similar to that observed 
for species like Er or Dy. This happens in spite of the 
fact that the MQDT model makes no assumptions about 
random couplings as in Random Matrix Theory. The 
differences between the predictions of our model and the 
Wigner-Dyson statistics are best seen when analyzing the 
number variance, which for our case is well described by 
a semi-Poisson distribution, even for very strongly over¬ 
lapping resonances. Thus, the current observational and 
theoretical evidence on the level spacing and its num¬ 
ber variance [30] combined with our result suggest that a 
semi-Poisson distribution may provide a better approx¬ 
imation to these quantities than a Wigner-Dyson one. 
This is a question that could be explored in future re¬ 
search. 

Our results show that correlations between resonance 
positions in ultracold collisions can have a different source 
than a chaotic distribution of bound levels. The effect of 
overlapping of resonances itself introduces level repulsion 
which will affect the analysis of experimental data. An¬ 
other particular experimental and theoretical challenge 
lies in measuring and interpreting the resonance widths. 
Recent analysis of Er data [33], based on the width of 
three-body loss features, showed discrepancies with RMT 
predictions. However, the three-body recombination pro¬ 
cess in the vicinity of narrow overlapping resonances is 
very difficult to describe precisely m and it is not obvi¬ 
ous what information can be extracted from this about 
the resonance widths. 

The MQDT model presented here can be extended to 
take into account more detailed assumptions about the 
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interaction potentials and provide an analysis of realistic 
species complementary to the one provided by ab initio 
calculations. For example, ref. [T5] shows that special 
non-chaotic eigenstates can emerge from the dense set of 
chaotic states and play a dominant role in determining 
the near-threshold bound state and scattering proper¬ 
ties. Our MQDT method readily explains the effect of 
such unique non-chaotic states on the bound state and 
scattering properties within an overall dense set of over¬ 


lapping chaotic resonances. Overall, a variety of theoret¬ 
ical approaches is needed to fully understand the role of 
quantum chaos in ultracold collisions 
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